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Background. In order to reduce the side-effects of chemotherapy, combined chemotherapy-electroporation (elec- 
trochemotherapy) has been suggested. Electroporation, application of appropriate electric pulses to biological cells, 
can significantly enhance molecular uptake of cells due to formation of transient pores in the cell membrane. It was 
experimentally demonstrated that the efficiency of electroporation is under the control of electric pulse parameters. 
However, the theoretical basis for these experimental results is not fully explained. In order to predict the outcome of 
experiments and optimize the efficiency of electroporation before each treatment, we developed a model to inves- 
tigate the effect of pulse shape on efficiency of electroporation. 

Results. Our model is based on a developed chemical-kinetics scheme and trapezium barrier model, while self- 
consistency was taken into account. This model is further supplemented with a molecular transport model to acquire 
the molecular uptake of cells. The investigated pulse shapes in this study were unipolar rectangular pulses with dif- 
ferent rise and fall times, triangular, sinusoidal and bipolar rectangular pulses and also sinusoidal modulated unipolar 
pulses with different percentages of modulation. The obtained results from our modelling and simulations are in good 
agreement with previously published experimental results. 

Conclusions. We therefore conclude that this model can be used to predict the effects of arbitrarily shaped electro- 
poration pulses on cell membrane conductivity and molecular transport across the cell membrane. 
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Introduction 

Cancer is a leading cause of death around the 
world and for this reason it has drawn the at- 
tention of many researchers. Chemotherapy has 
been used for many years and is one of the most 
common treatments for cancer. Cytotoxic chemo- 
therapeutic drugs are usually hydrophilic with 
very low transport through the cell membrane 
and thus high doses of these drugs are needed for 
treatment. Therefore, while chemotherapy can be 
quite effective in treating certain cancers by inter- 
fering with the cancer cell's ability to grow or re- 
produce, chemotherapeutic drugs reach all parts 
of the body, not just the cancer cells. Because of 
this, there may be many side-effects such as nau- 
sea, blood cell deficiency, fatigue and loss of hair 
during treatment. 
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Two decades ago, electrochemotherapy was 
suggested for its use in clinical treatment of local- 
ized tumors. 2 Electrochemotherapy consists of elec- 
troporation and chemotherapy. Electroporation is 
a technique in which permeability of the plasma 
membrane increases transiently and reversibly 
with appropriate pulse parameters 3 " 6 and is nowa- 
days widely used not only in electrochemothera- 
py, 7 " 9 but also in biotechnology 10 ' 11 and in medical 
applications such as gene electrotransfer 12 ' 13 and 
transdermal drug delivery. 14 Electroporation in 
combination with chemotherapy can increase drug 
delivery into the cells and consequently drug doses 
and thus side-effects of chemotherapy can be re- 
duced. 10 ' 11 ' 15 ' 16 

It was experimentally demonstrated in a number 
of studies that the efficiency of electroporation 
is under the control of electric pulse parameters 
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such as pulse amplitude, duration, and shape. 
Optimization of electric field parameters for suc- 
cessful electroporation requires time-consuming 
and costly experiments for different experimental 
criteria unless an appropriate model for this phe- 
nomenon can be suggested. Although the known 
models were proposed in previous studies 23 " 28 , they 
are still unable to explain the effect of some param- 
eters such as pulse shape, pulse repetition frequen- 
cy and number of pulses on the molecular uptake 
enhancement of the cells under exposure. 

The suggested mechanism for electroporation 
consists of structural changes resulting in forma- 
tion of transient aqueous pores in the cell mem- 
brane. In order to reveal the exact mechanisms and 
dynamics of pore formation and closure and more 
importantly resealing of the membrane, theoretical 
models have drawn a great deal of attention. 

In our present study, we investigated the effect 
of pulse shape on the efficiency of electroporation 



and, therefore, electrochemotherapy using mod- 
elling and simulation. Our model was based on a 
chemical-kinetics scheme with two types of pores 29 
which have been recently confirmed. 30 We used de- 
veloped equations with field-dependent rate coef- 
ficients in order to obtain the pore distribution on 
the membrane. Besides, the conductivity of pores 
was defined based on a trapezium barrier model 
for the image forces. 31 A self-consistent set of equa- 
tions was used to consider all simultaneous chang- 
es. This model was supplemented with a molecu- 
lar transport model for a single cell to acquire the 
molecular uptake of cells. The investigated pulse 
shapes in this study were unipolar rectangular 
pulses with different rise and fall times, triangular, 
sinusoidal and bipolar rectangular pulses and also 
sinusoidal modulated unipolar pulses with dif- 
ferent percentages of modulation - all previously 
used in experimental studies. 



Material and methods 

Model description 

Modified chemical-kinetics model for electroporation 

When a cell is exposed to an external electric field, the induced transmembrane voltage (ITV) starts to 
increase based on the Laplace equation which leads to structural changes of the cell membrane. Based on 
a previously suggested 29 , and recently confirmed 30 kinetic model, in the first step the intact closed lipids 
(C) transform to tilted lipid headgroups (CI). In the second step, the prepores (PI) are formed and finally, 
in the last step the final pores (P2) are formed. The sequential reaction can be described by: 

c 4= cl + pl + p2 m 

The permeability of the PI state is negligibly small and P2 is predominantly responsible for molecular 
uptake. Pore formation and closure are denoted by k { and k A (i=l,2,3) rate coefficients, respectively. For 
simplicity, the rate coefficients k l7 k 2 and k 3 are considered equal (k 1 =k 2 =k 3 =k p ). 29 The governed rate laws of 
constituting steps for the scheme [1] are: 



^^ = - kp [C{rJ)] + k_\C\{r,t)\ 
at 

d[C1( f' t)] = -k ([Cl(r,t)]- [C(r,f)])-*-i [Cl(r,t)] + k_ 2 [Pl(r,t)] 
at 

d[P Y ,t)] = -k p i\P\{7, f )] - [Cl(r, t)]) - k_ 2 [Pl(r, 0] + k_ 3 [P2(f, t)] 



^fM = k p [Pl(r,t)]-k_ 3 [P2(r,t)} [2] 
at 

where t and r denote time and position, respectively. [C], [CI], [PI] and [P2] show normalized distribution 
of each membrane lipid state relative to the initial value of the closed state [C(r,0)]. 29 

Regarding the Van't Hoff relationship in electro-thermodynamics, the rate coefficient of pore formation 
can be obtained from: 29,31 
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where ITV is the potential difference between the outer and inner layer of the membrane, p is the mean 
volume change due to pore formation, £ 0 is the permittivity of the vacuum and £ w and £ L are dielectric 
constants of water and lipids, respectively. k B is the Boltzmann constant, d m is the thickness of the mem- 
brane and T is temperature. While the pore formation rate coefficient k p is electric field-dependent, the 
closure rate coefficients (k. y k 2 and k 3 ) are constant and independent of electric field strength. 29 

Whenever electroporation occurs, an increase in conductivity during the pulse is observed 32 which can 
be explained by the formation of pores in the cell membrane. Based on the trapezium barrier model for 
the image forces, the intrinsic pore conductivities o pi (i=l and 2 represents PI and P2 pores, respectively) 
are expressed as follows: 31,33 



[4] 



where 



o ov + <7 ir 



-exp 



RT 



and a nA = 1 ■ 



RT 



[5] 



In the above equations, o ex and o in are the extracellular and intracellular conductivities respectively, n 
is the geometrical parameter of the trapezium model for energy barrier, F is the Faraday constant and tpf m . 
is the intrinsic pore barrier potential. 

Therefore, conductivity of the membrane (a m ) can be obtained by: 

a m (r , 0 = a m0 + [Pl(r, t)] x a pA + [P2(r, t)] x a pa [6] 



where <7 m0 is the physiological/baseline conductivity of the membrane. Thus conductivity at each point on 
the membrane changes with time during and after the pulse, depending on pore distribution variations 
which affect ITV and in turn the distribution of pores. 

Transmembrane molecular transport model 

Based on previous studies 34 " 36 , we defined two distinct phases for the electroporated membrane and two 
related transport mechanisms: the first one is the porated phase [P2] with relatively fast relaxation due to 
pore closure according to Eq. [1]. The second phase is the memory phase [M] due to enhanced membrane 
perturbation and ruffling with quite slow relaxation 34 " 36 which returns to its baseline value with a dual 
exponential decay function: 29 

[M] = [P2] e {B exp(-^0 + (1 - B) exp(-^0) [7] 



where [P2] e is the normalized distribution of [P2] pores at the end of the pulse, k f and k s are decay rate 
coefficients for this second phase and B is a constant. 

The considered transport mechanisms for these two phases were interactive diffusion through the 
pores and endocytotic-like transport through the permeabilized area of the membrane. 34 " 37 Thus, the per- 
meability of the membrane can be written as the sum of two distinct contributions: 

P m (r,t) = {[P2(r,t)]D p /d m )+&M(r,tj\P r ld n ) [8] 



where D p and D r are the attributed diffusion coefficients for interactive transport and endocytotic-like 
transport, respectively. 

While the membrane is being permeabilized due to the electric field, the molecules pass through the 
membrane due to a concentration gradient. A quantitative description of diffusion is contained in Fick's 
first law. The total flux can be approximated by j = P m (c out - c in ), where c out and c in are the outer and inner 
concentrations adjacent to the membrane. The total number of molecules transported through the mem- 
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TABLE 1. Values of parameters used in simulations 



Parameter 


Symbol 


Value 


Ref. 


Membrane thickness 


dm 


5e-9 m 


29 


Extracellular conductivity 


<J ex 


0.14 S/m a 


29, 38 


Intracellular conductivity 


°in 


0.3 S/m b 


29, 39 


Initial conductivity of membrane 




5e-7 S/m 


40 


Extracellular permittivity 


£ 0 


7.1e-10 As/Vm 


29 


Intracellular permittivity 




7.1e-10 As/Vm 


29 


Membrane permittivity 




4.4e-ll As/Vm c 


29 


Water relative dielectric constant 




80 As/Vm 


29 


Lipid relative dielectric constant 


£i 


2 As/Vm 


29 


Free diffusion coefficient 


D 0 


5e-10 m 2 /s 


29 


Zero-field equilibrium constant 


Ko 


2e-2 


29 


Mean average aqueous pore volume 


AV p 


9e-27 m 3 


29 


Intrinsic barrier potential of PI state 


J) 
<Pim\ 


0.13 V 


31 


Intrinsic barrier potential of P2 state 


<Pim2 


0.084 V 


31 


A geometrical parameter 


n 


0.12 


33 


Decay rate coefficient for CI 


k, 


10 5 s" 1 


41-43 


Decay rate coefficient for PI pores 


k- 2 


2000 s" 1 


41-43 


Decay rate coefficient for P2 pores 


k-3 


2 s" 1 


41 


Decay rate coefficients for endocytotic-like process 


k f ,k s 


0.044, 0.003 s" 1 


29 



a This is for SMEM. The range of extracellular medium is quite large. b Reported between (0.2-0.55) S/m c Reported between (4.4-5)*10~ n As/Vm 



brane (N) was computed with integration of transported molecules through the cell membrane over time 
and the cell surface: 

N = N A j \ jdSdt [9] 

t=os 

where S is the surface of the cell membrane, T is the time at which the quantity of transported molecules 
is to be determined and N A is Avogadro's number. 



Construction of the model 

The simulations in this study were performed 
using the COMSOL 3.3 package (COMSOL Inc., 
Burlington, MA) based on the finite element meth- 
od. To construct the geometrical model, a spheri- 
cal cell with radius of 5.6 um was located between 
two virtual electrodes. Since incorporating an ex- 
tremely thin membrane is problematic in meshing 
and solving the problem, we assigned the bound- 
ary condition to the membrane. 38 We neglected the 
resting transmembrane voltage. The initial intra- 
cellular and extracellular concentrations of probe 
were set to 0 and 10 mM, respectively. The diffu- 
sion coefficients for interactive diffusion and for an 



induced endocytotic-like process are considered 
as D 0 /5 and D 0 /10000. These two values, however, 
depend considerably on the type and size of the 
transported molecules. The necessary parameters 
used in our simulations are given in Table 1. Our 
simulation was designed to solve the Laplace equa- 
tion considering all related equations in this model 
(Eq. [2], [3], [6]) taking into account self-consistency 
of parameters to find the distribution of pores on 
the cell membrane, spatially and temporally, and 
all related parameters such as ITV, cell membrane 
conductivity and permeability. Afterwards, the up- 
take of the cells for each different pulse shape was 
obtained. All simulations were performed on a PC 
(2.8 GHz Pentium IV processor, 3 GB RAM) and 
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FIGURE 1. The investigated rectangular pulse shapes in this FIGURE 2. The (A) bipolar rectangular, (B) triangular and (C) 

study with rise and fall times of (A) 2, (B) 10 and (C) 100 |js, sinusoidal pulses considered in this study. 

respectively. 




FIGURE 3. The sine-modulated 50 kHz unipolar pulses with (A) 10% and (B) 90% modulation investigated in this paper. 
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FIGURE 4. Temporal evolution of the overall membrane conductivity during the pulse for (A) unipolar and bipolar rectangular, triangular and 
sinusoidal pulses and (B) 10% and 90% sine-modulated unipolar pulses. 



each simulation lasted 3-25 minutes depending on 
the considered pulse shape and number of pulses 
in each train of pulse. 

The investigated pulse shapes in this study were 
unipolar rectangular pulses with different rise and 
fall times of 2, 10 and 100 us (Figure 1); triangular, 
sinusoidal and bipolar rectangular pulses (Figure 
2); and also sinusoidal modulated unipolar pulses 
with different percentages of modulation of 10% 
and 90% with 50 kHz frequency (Figure 3). 

Results 

Immediately after the smoothed step pulse is 
switched on, ITV starts to increase based on the 
Laplace equation and causes membrane structur- 
al changes initiation, which in turn results in the 
membrane conductivity increase according to Eq. 
[6]. The temporal behaviour of average conductivi- 
ties over the cell membrane due to application of 
the considered pulse shapes (Figures 1, 2 and 3) are 
shown in Figure 4. All pulses were considered to 
have a peak of 1 kV/cm and total duration of 1 ms. 



It can be observed in Figure 4 that the overall 
conductivity changes for unipolar and bipolar 
pulses have negligible differences. The reason for 
this fact is a very quick switch between positive 
and negative voltage, as well as ignoring the rest- 
ing voltage in this model. Besides, a comparison 
between conductivity increases due to rectangular, 
triangular and sinusoidal pulses was performed. 
Figure 4A shows that the largest and smallest 
changes were due to rectangular and triangular 
pulse shapes, respectively. Figure 4B shows that 
a conductivity change due to 10% modulation is 
higher than the 90% one but both are still lower 
than the rectangular pulse. 

The temporal behaviour of averaged cell mem- 
brane permeability for pulses in Figures 1 to 3 is 
illustrated in Figure 5. Permeability changes occur 
slowly. Therefore, for bipolar pulses and modu- 
lated pulses in which the fall and rise is very fast, 
there is not enough time for resealing of permeabil- 
ity which causes different behaviour for membrane 
permeability related to membrane conductivity. 
Based on Figure 5, we expect the order of efficien- 
cy of pulses of the same peak as follows: unipolar 




FIGURE 5. Temporal evolution of the overall membrane permeability during the pulse for (A) unipolar and bipolar rectangular, triangular and 
sinusoidal pulses and (B) 10% and 90% sine-modulated unipolar pulses. 
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Figure 6A shows the results of simulation for 
8 pulses of bipolar rectangular, sine and triangu- 
lar pulses. It shows that the rectangular pulses are 
more efficient than sinusoidal pulses which in turn 
are more efficient than triangular pulses. These re- 
sults are in good agreement with experimentally 
obtained results. 9 

As can be seen in Figure 2A for the bipolar 
pulses, the pulse switch from positive to negative 
takes place very fast. During the switch time, the 
pore creation rate and, therefore, membrane con- 
ductivity decrease. But due to very short time of 
switching related to pulse duration, these changes 
are negligible in comparison to the conductivity 
change related to the unipolar pulse (Figure 4A). 
Consequently, the uptake due to unipolar pulses 
is larger than bipolar pulses but this difference is 
negligible and not observable (data not shown). 
While our simulation shows no significant differ- 
ence between these two pulse types, in experimen- 
tal results bipolar pulses are significantly more 
efficient than unipolar pulses. The reason for this 
inconsistency is most probably due to neglecting 
resting voltage in the simulations. 

In addition, Figure 6B demonstrates the com- 
parison between unipolar pulses of 0, 10 and 90% 
modulation. The results are also in good qualita- 
tive agreement with previously obtained experi- 
mental results. 9 The uptake enhancement results 
for 8 pulses of unipolar trapezoidal pulses of 1 ms 
duration with 2, 10 and 100 us rise and fall times 
are shown in Figure 6C. It can be seen from the fig- 
ure that there is no significant difference between 
these pulses which is again in good agreement 
with previously published experimental results. 19 



FIGURE 6. Dye uptake 16 minutes after pulse cessation versus electric field 
amplitude for (A) 8 pulses of 1 ms and 1 Hz unipolar and bipolar rectangu- 
lar, triangular and sinusoidal pulses, (B) one 1 ms pulse of 0%, 10% and 90% 
sine-modulated and (C) 8 pulses of 1 ms rectangular pulses with rise and 
fall times of 2, 10 and 100 ms. 



and bipolar rectangular, 10% modulated, sine, 90% 
modulated and finally bipolar triangular pulses. 
Note that unlike membrane conductivity in Figure 
4, membrane permeability does not recover as fast 
after the pulse ceases. 

To be able to check the validity of our simulation 
results, the uptake enhancement of the cell was 
calculated for the same pulse parameters of previ- 
ously obtained experimental results. 9 The chosen 
parameters were 8 pulses of 100 us duration and 1 
Hz pulse repetition frequency with different pulse 
strengths for each pulse shape. 



Conclusions 

The described model enables determination and 
prediction of all electrical and diffusion parameters 
for different pulse shapes. Thus, knowing electrical 
and diffusion properties of the cells and the specific 
dye, optimization of the electroporation protocol 
can be performed before the treatment. Our results 
show that rectangular pulses are more effective 
than the sinusoidal and triangular pulses. Besides, 
our results indicate that the higher the percentage of 
unipolar pulses modulation with sine shape pulses 
of 50 kHz, the lower the uptake enhancement of the 
cells. Moreover, the rise and fall times of unipolar 
rectangular pulses do not significantly affect the 
uptake of molecules by the cells. Our simulation re- 
sults are consistent with experimental observations. 
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